package grid_operation;

import org.gdal.gdal.Dataset;
import org.gdal.gdal.WarpOptions;
import org.gdal.gdal.gdal;
import org.gdal.gdalconst.gdalconstConstants;
import org.gdal.ogr.*;
import utils.SelectFile;

import java.util.ArrayList;
import java.util.List;
import java.util.Vector;

public class LimitedScopeGrid {
    public static void main(String[] args) {
        String tifTuCenMax = "E:\\TestInput\\TestInput2\\tucenghoudu.tif";
        String tifLanUse = "E:\\TestInput\\TestInput2\\landuse.tif";
        String outPath = "E:\\Test\\result\\";
        String maskShape3 = "E:\\TestInput\\TestInput2\\shape\\test.shp";

        //多边形数组
        double[][] pointsArray = new double[][]{
                new double[]{
                        48.58734,  //lat
                        87.72319   //lng
                },
                new double[]{
                        48.01728,
                        87.28161
                },
                new double[]{
                        47.65378,
                        88.28755
                },
                new double[]{
                        48.39410,
                        88.86063
                },
                new double[]{
                        48.58734,
                        87.72319
                },
        };

        List<String> tifs = new ArrayList<>();
        tifs.add(tifTuCenMax);
        tifs.add(tifLanUse);

        List<String> files = limiteGrid(tifs, outPath, maskShape3, pointsArray);
    }

    /**
     * 切割tif算法
     * @param tifs 传入的tif集合
     * @param outPath 切割完成后输出tif的保存目录
     * @param maskShape 多边形生产的临时shape路径
     * @param pointsArray 多边形顶点数组
     * @return
     */
    private static List<String> limiteGrid(List<String> tifs, String outPath, String maskShape, double[][] pointsArray) {
        gdal.AllRegister();
        gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES");
        gdal.SetConfigOption("SHAPE_ENCODING", "UTF-8");

        if (tifs.isEmpty()) {
            return new ArrayList<>();
        }

        org.gdal.ogr.Driver driver = SelectFile.driverInit();
        Dataset hDataset = gdal.Open(tifs.get(0));
        DataSource oDS = driver.CreateDataSource(maskShape, null);
        Layer oLayer = oDS.CreateLayer("TestPolygon", hDataset.GetSpatialRef(), ogr.wkbPolygon, null);

        FieldDefn oFieldID = new FieldDefn("FieldID", ogr.OFTInteger);
        oLayer.CreateField(oFieldID);

        FeatureDefn oDefn = oLayer.GetLayerDefn();

        StringBuffer sb = new StringBuffer();
        int index = 0;
        for (double[] points : pointsArray) {
            if (index == 0) {
                sb.append("POLYGON ((");
                sb.append(points[1]);
                sb.append(" ");
                sb.append(points[0]);
            } else if (index == pointsArray.length - 1) {
                sb.append(",");
                sb.append(points[1]);
                sb.append(" ");
                sb.append(points[0]);
                sb.append("))");
            } else {
                sb.append(",");
                sb.append(points[1]);
                sb.append(" ");
                sb.append(points[0]);
            }
            index++;
        }

        // 创建多边形要素
        Feature oFeatureTriangle = new Feature(oDefn);
        oFeatureTriangle.SetField(0, 0);
        Geometry geomTriangle = Geometry.CreateFromWkt(sb.toString());
        oFeatureTriangle.SetGeometry(geomTriangle);
        oLayer.CreateFeature(oFeatureTriangle);

        oDS.SyncToDisk();

        List<String> clips = new ArrayList<>();
        for (String tif: tifs) {
            String outPathFile = outPath + "\\" + System.currentTimeMillis() + ".tif";
            clips.add(clip(tif,
                    outPathFile,
                    maskShape));
        }

        return clips;
    }


    /**
     * 执行切割
     * @param srcFile 原图tif
     * @param dstFile 输出tif
     * @param shpFile 切割shape
     * @return
     */
    public static String clip(String srcFile, String dstFile, String shpFile) {
        // 读取要切的原图
        Dataset srcDs = gdal.Open(srcFile, gdalconstConstants.GA_ReadOnly);
        if (srcDs == null) {
            System.err.println("GDALOpen failed - " + gdal.GetLastErrorNo());
            System.err.println(gdal.GetLastErrorMsg());
            return dstFile;
        }

        // 构建warp参数
        Vector<String> options = new Vector<String>();
        options.add("-of");
        options.add("GTiff");
        options.add("-cutline");//按照矢量图像线裁剪
        options.add(shpFile);
        options.add("-crop_to_cutline"); //裁剪后图层范围使用矢量图层大小
        options.add("-dstnodata");
        options.add("255"); //nulldata的值

        WarpOptions wo = new WarpOptions(options);
        Dataset warp = gdal.Warp(dstFile, new Dataset[]{srcDs}, wo);

        // 获取切图的原点坐标信息
        double[] warpTransform = warp.GetGeoTransform();
        System.out.println("warp原点坐标   = " + warpTransform[0] + "," + warpTransform[3]);
        System.out.println("warp像素坐标差 = " + warpTransform[1] + "," + warpTransform[5]);

        srcDs.delete();

        return dstFile;
    }
}
